ValueFolderTemp='/sps/isc/rneveu/EtudeIRM/Nature/ResultsJournal/';
ValueFolderTemp2='/sps/isc/rneveu/EtudeIRM/Nature/ResultsEditor/';
load([ValueFolderTemp 'ResultsThirdLevelDistributionOutliers.mat']);
VecDistribution=VecResults;
load('/sps/isc/rneveu/EtudeIRM/Nature/EtudeNature.mat');
VecEditors1Temp(1,:)=[];
VecValueEditor=VecEditors1Temp;
load('/sps/isc/rneveu/EtudeIRM/Nature/DonneesNature.mat');
VecValueJournal=VecJournalTemp;
ValueCutoff=1.25;
VecResultsTemp1=[];
VecResultsTemp2=[];
VecResultsTemp13=[];
VecResultsTemp23=[];
VecResultsTemp11=[];
VecResultsTemp22=[];
VecValueYear=[0 1 2 3 4 5 10];
VecValueMax=[80 8 15 23 30 38 76];
nMonthMax=12;
n1Temp=0;
n1Temp3=0;
VecFigure=[];
VecFigure2=[];
VecValueFigure=[];
VecValueSlope=[1.34 1.98];
for n1=1:323
    ValueTemp=find(VecDistribution(:,1)==n1);
    ValueStartEditor=VecValueEditor{n1,22};
    if size(ValueTemp,1)==1
        if VecDistribution(ValueTemp,2)<=ValueCutoff
            ValueJournal=VecValueEditor{n1,11};
            n3=1;
            TestTemp=0;
            while n3<=size(VecValueJournal,1)&&TestTemp==0
                if strcmpi(ValueJournal,VecValueJournal{n3,1})==1
                    TestTemp=1;
                else
                    n3=n3+1;
                end
            end
            load([ValueFolderTemp2 'FileAnalysisEditor' num2str(n1) '.mat']);
            load([ValueFolderTemp2 'FileAnalysisCharacteristicsEditor' num2str(n1) '.mat']);
            VecValueTemp2=[];
            VecValueTemp22=[];
            n2Temp=0;
            if n3>=32
                ValueSlope=VecValueSlope(1,1);
            else
                ValueSlope=VecValueSlope(1,2);
            end
            for n2=1:33
                if n2==n3
%                     eval(['VecValueTemp=VecResultsEditorAuthors2FilteredJournal' num2str(n2) ';']);
%                     eval(['VecValueTemp=VecResultsEditorAuthors2Journal' num2str(n2) ';']);
                    eval(['VecValueTemp=VecResultsEditorAffiliation2Journal' num2str(n2) ';']);
                    if sum(VecValueTemp(2,1:24))>0
                        VecValueTemp22=VecValueTemp;
                    else
                        VecValueTemp22=NaN(2,168);
                    end
                else
%                     eval(['VecValueTemp=VecResultsEditorAuthors1FilteredJournal' num2str(n2) ';']);
%                     eval(['VecValueTemp=VecResultsEditorAuthors1Journal' num2str(n2) ';']);
                    eval(['VecValueTemp=VecResultsEditorAffiliation1Journal' num2str(n2) ';']);
                    if sum(VecValueTemp(2,1:24))>0
                        VecValueTemp2=[VecValueTemp2;VecValueTemp];
                        n2Temp=n2Temp+1;
                    end
                end
            end
            VecValueTemp23=1:2:n2Temp*2;
            VecValueTemp23=VecValueTemp23';
            for n4=1:size(VecValueTemp23,1)
                for n5=2:2:size(VecValueTemp2,2)/2
                    if VecValueTemp2(VecValueTemp23(n4,1)+1,(n5-1)*2+1)==0||VecValueTemp2(VecValueTemp23(n4,1)+1,(n5-1)*2+2)==0
                        VecValueTemp2(VecValueTemp23(n4,1),(n5-1)*2+1)=NaN;
                        VecValueTemp2(VecValueTemp23(n4,1)+1,(n5-1)*2+1)=NaN;
                        VecValueTemp2(VecValueTemp23(n4,1),(n5-1)*2+2)=NaN;
                        VecValueTemp2(VecValueTemp23(n4,1)+1,(n5-1)*2+2)=NaN;
                    end
                end
            end
            for n5=2:2:size(VecValueTemp22,2)/2
                if VecValueTemp22(2,(n5-1)*2+1)==0||VecValueTemp22(2,(n5-1)*2+2)==0
                    VecValueTemp22(1,(n5-1)*2+1)=NaN;
                    VecValueTemp22(2,(n5-1)*2+1)=NaN;
                    VecValueTemp22(1,(n5-1)*2+2)=NaN;
                    VecValueTemp22(2,(n5-1)*2+2)=NaN;
                end
            end
            n5=2*nMonthMax+8;
            VecValueTemp24=2:2:n2Temp*2;
            VecValueTemp24=VecValueTemp24';
%             VecValueFigure=[VecValueFigure;(VecValueTemp2(VecValueTemp23,(n5-1)*2+1)-VecValueTemp2(VecValueTemp23,(n5-1)*2+2))*2./(VecValueTemp2(VecValueTemp23,(n5-1)*2+1)+VecValueTemp2(VecValueTemp23,(n5-1)*2+2))*100 VecValueTemp2(VecValueTemp23,(n5-1)*2+1)./VecValueTemp2(VecValueTemp24,(n5-1)*2+1)*100-VecValueTemp2(VecValueTemp23,(n5-1)*2+2)./VecValueTemp2(VecValueTemp24,(n5-1)*2+2)*100 VecResultsEditor(1,3)*ones(size(VecValueTemp23,1),1) VecResultsEditor(1,27)*ones(size(VecValueTemp23,1),1) n1*ones(size(VecValueTemp23,1),1) ones(size(VecValueTemp23,1),1)];
%             VecValueFigure=[VecValueFigure;(VecValueTemp22(1,(n5-1)*2+1)-VecValueTemp22(1,(n5-1)*2+2))*2./(VecValueTemp22(1,(n5-1)*2+1)+VecValueTemp22(1,(n5-1)*2+2))*100 VecValueTemp22(1,(n5-1)*2+1)./VecValueTemp22(2,(n5-1)*2+1)*100-VecValueTemp22(1,(n5-1)*2+2)./VecValueTemp22(2,(n5-1)*2+2)*100 VecResultsEditor(1,3) VecResultsEditor(1,27) n1 2];
            VecValueTemp25=VecValueTemp2(VecValueTemp23,:)./VecValueTemp2(VecValueTemp24,:);
            VecValueTemp25=nanmean(VecValueTemp25,1);
            VecValueTemp26=nanmean(VecValueTemp2(VecValueTemp23,:)-VecValueTemp2(VecValueTemp24,:),1);
            VecValueTemp27=nanmean(VecValueTemp2(VecValueTemp23,:),1);
            if size(VecValueTemp27,2)>0
                % Case of articles of co authors which are published in all Nature journals except the editor's one
                VecValueTemp3=[];
                for n2=1:size(VecValueYear,2)*nMonthMax
%                     VecValueTemp3=[VecValueTemp3;(VecValueTemp27(1,(n2-1)*2+1)-VecValueTemp27(1,(n2-1)*2+2))*2/(VecValueTemp27(1,(n2-1)*2+1)+VecValueTemp27(1,(n2-1)*2+2))*100 VecValueTemp22(1,(n2-1)*2+1)/VecValueTemp22(2,(n2-1)*2+1)*100-(VecValueTemp22(1,(n2-1)*2+2)/VecValueTemp22(2,(n2-1)*2+2))*100-VecValueTemp25(1,(n2-1)*2+1)*100+VecValueTemp25(1,(n2-1)*2+2)*100];
                    VecValueTemp3=[VecValueTemp3;VecValueTemp27(1,(n2-1)*2+1)-VecValueTemp27(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1)-VecValueTemp22(1,(n2-1)*2+2)];
                    VecResultsTemp11=[VecResultsTemp11;VecValueTemp27(1,(n2-1)*2+1) VecValueTemp27(1,(n2-1)*2+2) 1 VecResultsEditor(1,23) VecResultsEditor(1,3) VecResultsEditor(1,16) VecResultsEditor(1,27)];
                end
                VecValueTemp4=[VecResultsEditor(1,27) VecResultsEditor(1,1)+VecResultsEditor(1,2) VecResultsEditor(1,16) VecResultsEditor(1,27)];
                VecValueTemp3=[VecValueTemp3 repmat(VecValueTemp4,size(VecValueTemp3,1),1) n1*ones(size(VecValueTemp3,1),1)];
                if size(ValueStartEditor,2)<3
                    VecResultsTemp1=[VecResultsTemp1;VecValueTemp3];
                else
                    VecResultsTemp13=[VecResultsTemp13;VecValueTemp3];
                end
                % Case of articles of co authors which are published in the journal of the editor only
                VecValueTemp32=[];
                for n2=1:size(VecValueYear,2)*nMonthMax
%                     VecValueTemp32=[VecValueTemp32;(VecValueTemp22(1,(n2-1)*2+1)-VecValueTemp22(1,(n2-1)*2+2))*2/(VecValueTemp22(1,(n2-1)*2+1)+VecValueTemp22(1,(n2-1)*2+2))*100 VecValueTemp22(1,(n2-1)*2+1)/VecValueTemp22(2,(n2-1)*2+1)*100-(VecValueTemp22(1,(n2-1)*2+2)/VecValueTemp22(2,(n2-1)*2+2))*100];
                    VecValueTemp32=[VecValueTemp32;VecValueTemp22(1,(n2-1)*2+1)-VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1)./VecValueTemp22(2,(n2-1)*2+1)*100-VecValueTemp22(1,(n2-1)*2+2)./VecValueTemp22(2,(n2-1)*2+2)*100];
                    VecResultsTemp22=[VecResultsTemp22;VecValueTemp22(1,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) 1 VecResultsEditor(1,23) VecResultsEditor(1,3) VecResultsEditor(1,16) VecResultsEditor(1,27)];
                end
                VecValueTemp32=[VecValueTemp32 repmat(VecValueTemp4,size(VecValueTemp3,1),1) n1*ones(size(VecValueTemp32,1),1)];
                if size(ValueStartEditor,2)<3
                    VecResultsTemp2=[VecResultsTemp2;VecValueTemp32];
                    n2=2*nMonthMax+8;
                    VecFigure=[VecFigure;VecValueTemp27(1,(n2-1)*2+2) VecValueTemp27(1,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,3) VecResultsEditor(1,27) n1 1 1];
                    VecFigure2=[VecFigure2;VecValueTemp22(2,(n2-1)*2+2) VecValueTemp22(2,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,1) VecResultsEditor(1,16) VecResultsEditor(1,27) ValueSlope n1 n3 1 1];
                else
                    VecResultsTemp23=[VecResultsTemp23;VecValueTemp32];
                    n2=2*nMonthMax+8;
                    VecFigure=[VecFigure;VecValueTemp27(1,(n2-1)*2+2) VecValueTemp27(1,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,3) VecResultsEditor(1,27) n1 1 2];
                    VecFigure2=[VecFigure2;VecValueTemp22(2,(n2-1)*2+2) VecValueTemp22(2,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,1) VecResultsEditor(1,16) VecResultsEditor(1,27) ValueSlope n1 n3 1 2];
                end
                VecValueTemp2=[];
                VecValueTemp22=[];
                n2Temp=0;
                for n2=1:33
                    if n2==n3
                        eval(['VecValueTemp=VecResultsEditorAuthors2FilteredJournal' num2str(n2) ';']);
%                         eval(['VecValueTemp=VecResultsEditorAffiliation2Journal' num2str(n2) ';']);
                        VecValueTemp22=VecValueTemp;
                    else
                        eval(['VecValueTemp=VecResultsEditorAuthors1FilteredJournal' num2str(n2) ';']);
%                         eval(['VecValueTemp=VecResultsEditorAffiliation1Journal' num2str(n2) ';']);
                        VecValueTemp2=[VecValueTemp2;VecValueTemp];
                        n2Temp=n2Temp+1;
                    end
                end
                VecValueTemp23=1:2:n2Temp*2;
                VecValueTemp23=VecValueTemp23';
                for n4=1:size(VecValueTemp23,1)
                    for n5=1:size(VecValueTemp2,2)
                        if VecValueTemp2(VecValueTemp23(n4,1)+1,n5)==0
                            VecValueTemp2(VecValueTemp23(n4,1),n5)=NaN;
                            VecValueTemp2(VecValueTemp23(n4,1)+1,n5)=NaN;
                        end
                    end
                end
                for n5=1:size(VecValueTemp22,2)
                    if VecValueTemp22(2,n5)==0
                        VecValueTemp22(1,n5)=NaN;
                        VecValueTemp22(2,n5)=NaN;
                    end
                end
                VecValueTemp24=2:2:n2Temp*2;
                VecValueTemp24=VecValueTemp24';
                VecValueTemp25=nanmean(VecValueTemp25,1);
                VecValueTemp26=nanmean(VecValueTemp2(VecValueTemp23,:)-VecValueTemp2(VecValueTemp24,:),1);
                VecValueTemp27=nanmean(VecValueTemp2(VecValueTemp23,:),1);
                if size(ValueStartEditor,2)<3
                    n2=2*nMonthMax+8;
                    VecFigure=[VecFigure;VecValueTemp27(1,(n2-1)*2+2) VecValueTemp27(1,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,3) VecResultsEditor(1,27) n1 2 1];
                    VecFigure2=[VecFigure2;VecValueTemp22(2,(n2-1)*2+2) VecValueTemp22(2,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,1) VecResultsEditor(1,16) VecResultsEditor(1,27) ValueSlope n3 n1 2 1];
                else
                    n2=2*nMonthMax+8;
                    VecFigure=[VecFigure;VecValueTemp27(1,(n2-1)*2+2) VecValueTemp27(1,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,3) VecResultsEditor(1,27) n1 2 2];
                    VecFigure2=[VecFigure2;VecValueTemp22(2,(n2-1)*2+2) VecValueTemp22(2,(n2-1)*2+1) VecValueTemp22(1,(n2-1)*2+2) VecValueTemp22(1,(n2-1)*2+1) VecResultsEditor(1,1) VecResultsEditor(1,16) VecResultsEditor(1,27) ValueSlope n3 n1 2 2];
                end
                if size(ValueStartEditor,2)<3
                    n1Temp=n1Temp+1;
                else
                    n1Temp3=n1Temp3+1;
                end
            end
        end
    end
end

% save([ValueFolderTemp2 'FichierFigureCoAuthorsAffiliationOtherJournalsNonFiltered.mat'],'VecFigure');
% save([ValueFolderTemp2 'FichierFigureCoAuthorsAffiliationOtherJournals.mat'],'VecFigure');
save([ValueFolderTemp2 'FichierFigureCoAuthorsAffiliationOtherJournals2.mat'],'VecFigure2');
% save([ValueFolderTemp2 'AnalyseMixedModelCoAuthorsNatureNonFiltered.mat'],'VecValueFigure');
% save([ValueFolderTemp2 'AnalyseRobustnessCoAuthorsAffiliationNatureNonFiltered.mat'],'VecResultsTemp1','VecResultsTemp23','n1Temp3','VecResultsTemp2','VecValueYear','nMonthMax','n1Temp');


VecResults=[];
VecResults2=[];
VecResults1=[];
VecResults3=[];
VecResults4=[];
% n1Temp=n1Temp+n1Temp3;
% VecResultsTemp1=[VecResultsTemp1;VecResultsTemp13];
% VecResultsTemp2=[VecResultsTemp2;VecResultsTemp23];
for n1=1:size(VecValueYear,2)
    for n2=1:nMonthMax
        % Articles of co authors which are published in all Nature journals except the editor's one
        VecResultsTemp=[];
        for n3=1:n1Temp
            VecResultsTemp=[VecResultsTemp;VecResultsTemp1((n3-1)*size(VecValueYear,2)*nMonthMax+(n1-1)*nMonthMax+n2,:)];
        end
        ValueTemp=find(VecResultsTemp(:,1)<Inf&isnan(VecResultsTemp(:,1))==0);
        ValueTemp2=find(VecResultsTemp(:,2)<Inf&isnan(VecResultsTemp(:,2))==0);
        ValueTemp32=find(abs(VecResultsTemp(ValueTemp,1)-mean(VecResultsTemp(ValueTemp,1)))<=3*sqrt(var(VecResultsTemp(ValueTemp,1))));
        VecValueResultsTemp=VecResultsTemp;
        [ValueResults1Temp,ValueResults2]=signrank(VecResultsTemp(ValueTemp(ValueTemp32,1),1));
        [ValueResults3Temp,ValueResults4]=signrank(VecResultsTemp(ValueTemp(ValueTemp32,1),2));
        [ValueResults1,ValueResults2]=ttest(VecResultsTemp(ValueTemp(ValueTemp32,1),1));
        [ValueResults3,ValueResults4]=ttest(VecResultsTemp(ValueTemp(ValueTemp32,1),2));
%         [ValueResultsTemp,stats1]=robustfit(VecResultsTemp(ValueTemp,3),VecResultsTemp(ValueTemp,1));
%         [ValueResultsTemp2,stats2]=robustfit(VecResultsTemp(ValueTemp,4),VecResultsTemp(ValueTemp,1));
%         [ValueResultsTemp3,stats3]=robustfit(VecResultsTemp(ValueTemp,6),VecResultsTemp(ValueTemp,1));
%         [ValueResults2Temp,stats4]=robustfit(VecResultsTemp(ValueTemp2,3),VecResultsTemp(ValueTemp2,2));
%         [ValueResults2Temp2,stats5]=robustfit(VecResultsTemp(ValueTemp2,4),VecResultsTemp(ValueTemp2,2));
%         [ValueResults2Temp3,stats6]=robustfit(VecResultsTemp(ValueTemp2,6),VecResultsTemp(ValueTemp2,2));
        [ValueResultsTemp,ValueTemp3,stats1]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),1),{VecResultsTemp(ValueTemp(ValueTemp32,1),3)},'continuous',[1],'model','full','display','off');
        [ValueResultsTemp2,ValueTemp3,stats2]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),1),{VecResultsTemp(ValueTemp(ValueTemp32,1),4)},'continuous',[1],'model','full','display','off');
        [ValueResultsTemp3,ValueTemp3,stats3]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),1),{VecResultsTemp(ValueTemp(ValueTemp32,1),5)},'continuous',[1],'model','full','display','off');
        [ValueResults2Temp,ValueTemp3,stats4]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),2),{VecResultsTemp(ValueTemp(ValueTemp32,1),3)},'continuous',[1],'model','full','display','off');
        [ValueResults2Temp2,ValueTemp3,stats5]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),2),{VecResultsTemp(ValueTemp(ValueTemp32,1),4)},'continuous',[1],'model','full','display','off');
        [ValueResults2Temp3,ValueTemp3,stats6]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),2),{VecResultsTemp(ValueTemp(ValueTemp32,1),5)},'continuous',[1],'model','full','display','off');
        VecResults=[VecResults;mean(VecResultsTemp(ValueTemp(ValueTemp32,1),1)) sqrt(var(VecResultsTemp(ValueTemp(ValueTemp32,1),1))) ValueResults2 ValueResults1Temp mean(VecResultsTemp(ValueTemp(ValueTemp32,1),2)*100) sqrt(nanvar(VecResultsTemp(ValueTemp(ValueTemp32,1),2)*100)) ValueResults4 ValueResults3Temp n2 VecValueYear(1,n1) ValueResultsTemp' stats1.coeffs' NaN ValueResultsTemp2' stats2.coeffs' NaN ValueResultsTemp3' stats3.coeffs' NaN ValueResults2Temp' stats4.coeffs' NaN ValueResults2Temp2' stats5.coeffs' NaN ValueResults2Temp3' stats6.coeffs' NaN size(ValueTemp,1) size(ValueTemp2,1)];
        % Articles of co authors which are published in the editor's journal only
        VecResultsTemp=[];
        for n3=1:n1Temp
            VecResultsTemp=[VecResultsTemp;VecResultsTemp2((n3-1)*size(VecValueYear,2)*nMonthMax+(n1-1)*nMonthMax+n2,:)];
        end
        VecResultsTemp2Temp=[];
        for n3=1:n1Temp3
            VecResultsTemp2Temp=[VecResultsTemp2Temp;VecResultsTemp23((n3-1)*size(VecValueYear,2)*nMonthMax+(n1-1)*nMonthMax+n2,:)];
        end
        ValueTemp=find(VecResultsTemp(:,1)<Inf&isnan(VecResultsTemp(:,1))==0);
        ValueTemp2=find(VecResultsTemp(:,2)<Inf&isnan(VecResultsTemp(:,2))==0);
        ValueTemp32=find(abs(VecResultsTemp(ValueTemp,1)-mean(VecResultsTemp(ValueTemp,1)))<=3*sqrt(var(VecResultsTemp(ValueTemp,1))));
        ValueTemp22=find(VecResultsTemp2Temp(:,2)<Inf&isnan(VecResultsTemp2Temp(:,2))==0);
%         if n1==3&&n2==8
%             VecValueAffiliation=VecResultsTemp(:,1);
% %             VecValueCoAuthors=VecResultsTemp(:,1);
%         end
        [ValueResults1Temp,ValueResults2]=signrank(VecResultsTemp(ValueTemp(ValueTemp32,1),1));
        [ValueResults3Temp,ValueResults4]=signrank(VecResultsTemp(ValueTemp(ValueTemp32,1),2));
        [ValueResults1,ValueResults2]=ttest(VecResultsTemp(ValueTemp(ValueTemp32,1),1)); % Relative difference in the number of co authors' articles
        [ValueResults3,ValueResults4]=ttest(VecResultsTemp(ValueTemp(ValueTemp32,1),2)); % Difference in rates of co authors articles (N co authors' articles after/N total articles after -N co authors articles before/N total articles before)
        [ValueResults3Temp2,ValueResults4Temp]=ttest2(VecResultsTemp(ValueTemp2,1),VecResultsTemp2Temp(ValueTemp22,1));
%         [ValueResultsTemp,stats1]=robustfit(VecResultsTemp(ValueTemp,3),VecResultsTemp(ValueTemp,1));
%         [ValueResultsTemp2,stats2]=robustfit(VecResultsTemp(ValueTemp,4),VecResultsTemp(ValueTemp,1));
%         [ValueResultsTemp3,stats3]=robustfit(VecResultsTemp(ValueTemp,6),VecResultsTemp(ValueTemp,1));
%         [ValueResults2Temp,stats4]=robustfit(VecResultsTemp(ValueTemp2,3),VecResultsTemp(ValueTemp2,2));
%         [ValueResults2Temp2,stats5]=robustfit(VecResultsTemp(ValueTemp2,4),VecResultsTemp(ValueTemp2,2));
%         [ValueResults2Temp3,stats6]=robustfit(VecResultsTemp(ValueTemp2,6),VecResultsTemp(ValueTemp2,2));
        [ValueResultsTemp,ValueTemp3,stats1]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),1),{VecResultsTemp(ValueTemp(ValueTemp32,1),3)},'continuous',[1],'model','full','display','off');
        [ValueResultsTemp2,ValueTemp3,stats2]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),1),{VecResultsTemp(ValueTemp(ValueTemp32,1),4)},'continuous',[1],'model','full','display','off');
        [ValueResultsTemp3,ValueTemp3,stats3]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),1),{VecResultsTemp(ValueTemp(ValueTemp32,1),5)},'continuous',[1],'model','full','display','off');
        [ValueResults2Temp,ValueTemp3,stats4]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),2),{VecResultsTemp(ValueTemp(ValueTemp32,1),3)},'continuous',[1],'model','full','display','off');
        [ValueResults2Temp2,ValueTemp3,stats5]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),2),{VecResultsTemp(ValueTemp(ValueTemp32,1),4)},'continuous',[1],'model','full','display','off');
        [ValueResults2Temp3,ValueTemp3,stats6]=anovan(VecResultsTemp(ValueTemp(ValueTemp32,1),2),{VecResultsTemp(ValueTemp(ValueTemp32,1),5)},'continuous',[1],'model','full','display','off');
        VecResults2=[VecResults2;mean(VecResultsTemp(ValueTemp(ValueTemp32,1),1)) sqrt(nanvar(VecResultsTemp(ValueTemp(ValueTemp32,1),1))) ValueResults2 ValueResults1Temp mean(VecResultsTemp(ValueTemp(ValueTemp32,1),2)) sqrt(nanvar(VecResultsTemp(ValueTemp(ValueTemp32,1),2))) median(VecResultsTemp(ValueTemp(ValueTemp32,1),2)) quantile(VecResultsTemp(ValueTemp(ValueTemp32,1),2),0.25) quantile(VecResultsTemp(ValueTemp(ValueTemp32,1),2),0.75) ValueResults4 ValueResults3Temp n2 VecValueYear(1,n1) ValueResultsTemp' stats1.coeffs' NaN ValueResultsTemp2' stats2.coeffs' NaN ValueResultsTemp3' stats3.coeffs' NaN ValueResults2Temp' stats4.coeffs' NaN ValueResults2Temp2' stats5.coeffs' NaN ValueResults2Temp3' stats6.coeffs' NaN size(ValueTemp,1) size(ValueTemp2,1)];
        VecResults4=[VecResults4;mean(VecResultsTemp(ValueTemp2,2)) sqrt(nanvar(VecResultsTemp(ValueTemp2,2)))/sqrt(size(ValueTemp2,1)) mean(VecResultsTemp2Temp(ValueTemp22,2)) sqrt(nanvar(VecResultsTemp2Temp(ValueTemp22,2)))/sqrt(size(ValueTemp22,1)) ValueResults4Temp n2 VecValueYear(1,n1)];
        VecValueResultsTemp=VecResultsTemp(:,1)-VecValueResultsTemp(:,1);
        ValueTemp=find(VecValueResultsTemp(:,1)<Inf&isnan(VecValueResultsTemp(:,1))==0);
        [ValueResults1Temp,ValueResults2]=signrank(VecValueResultsTemp(ValueTemp,1));
        [ValueResults1,ValueResults2]=ttest(VecValueResultsTemp(ValueTemp,1));
        [ValueResults2Temp2,ValueTemp3,stats5]=anovan(VecValueResultsTemp(ValueTemp,1),{VecResultsTemp(ValueTemp,6)},'continuous',[1],'model','full','display','off');
        VecResults3=[VecResults3;nanmean(VecValueResultsTemp(ValueTemp,1)) sqrt(nanvar(VecValueResultsTemp(ValueTemp,1))) ValueResults1Temp ValueResults2 n2 VecValueYear(1,n1) ValueResults2Temp2' stats5.coeffs'];
        VecResultsTemp=[];
        for n3=1:n1Temp
%             VecResultsTemp=[VecResultsTemp;VecResultsTemp22(((n3-1)*size(VecValueYear,2)*nMonthMax+(n1-1)*nMonthMax+n2-1)*2+1:((n3-1)*size(VecValueYear,2)*nMonthMax+(n1-1)*nMonthMax+n2-1)*2+2,:)];
            VecResultsTemp=[VecResultsTemp;VecResultsTemp22((n3-1)*size(VecValueYear,2)*nMonthMax+(n1-1)*nMonthMax+n2,:)];
        end
%         ValueResults=anovan(VecResultsTemp(:,1),{VecResultsTemp(:,2) VecResultsTemp(:,3) VecResultsTemp(:,4)},'continuous',[1 3],'model',[1 0 0;0 1 0;0 0 1;1 1 0;1 0 1;0 1 1;1 1 1],'display','off');
        ValueResults=anovan(VecResultsTemp(:,1),{VecResultsTemp(:,2) VecResultsTemp(:,4)},'continuous',[1 2],'model',[1 0;0 1;1 1],'display','off');
        VecResults1=[VecResults1;ValueResults' n2 VecValueYear(1,n1)];
    end
end

% ValueTemp=find(VecResults(:,10)==0);
% VecResults(ValueTemp,10)=31*ones(size(ValueTemp,1),1);
% VecValueTemp=VecResults(:,10)*12+VecResults(:,9);
% [VecValueTemp,ValueTemp2]=sort(VecValueTemp);
% VecResults=VecResults(ValueTemp2,:);
% VecResults2=VecResults2(ValueTemp2,:);
% figure
% plot(VecValueTemp,VecResults2(:,5),'r-','LineWidth',4)
% hold on
% errorbar(VecValueTemp,VecResults2(:,5),VecResults2(:,6),'k+','LineWidth',1)
% hold on
% plot(VecValueTemp,VecResults(:,5),'-','LineWidth',4)
% hold on
% errorbar(VecValueTemp,VecResults(:,5),VecResults(:,6),'k+','LineWidth',1)



        
